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ABSTRACT 

We construct a physically motivated analytical model for the quasar luminosity func- 
tion, based on the joint star formation and feeding of massive black holes suggested 
by the observed correlation between the black hole mass and the stellar mass of the 
hosting spheroids. The parallel growth of massive black holes and host galaxies is 
assumed to be triggered by major mergers of haloes. The halo major merger rate is 
computed in the frame of the extended Press-Schechter model. The evolution of black 
holes on cosmological timescales is achieved by the integration of the governing set 
of differential equations, established from a few reasonable assumptions that account 
for the distinct (Eddington-limited or supply-limited) accretion regimes. Finally, the 
typical lightcurves of the reactivated quasars are obtained under the assumption that, 
in such accretion episodes, the fall of matter onto the black hole is achieved in a self- 
regulated stationary way. The predicted quasar luminosity function is compared to 
the luminosity functions of the 2dF QSO sample and other, higher redshift data. We 
find good agreement in all cases, except for z < 1 where the basic assumption of our 
model is likely to break down. 

Key words: cosmology - dark-matter - galaxies: formation, evolution - active galac- 
tic nuclei: radiation model, evolution - quasars: luminosity function, evolution 



1 INTRODUCTION 

In the last decade, a considerable progress has been achieved 
towards the understanding of the nature and evolution of 
Active Galactic Nuclei (AGN). It is now commonly accepted 
that they are powered by Massive Black Holes (MBHs) 
fuelled by material originating in the hosting dark-matter 
haloes and infalling from the close galactic environment in 
the process described as accretion. New models connecting 
the AGN phenomenon to the hierarchical growth of haloes 
have been developed, based on increasingly robust theoret- 
ical background substituting phenomenology. In this frame, 
the masses of MBHs are usually assumed to scale with the 
masses of the hosting dark-matter haloes with density es- 
timated from the Press & Schechter (1974) formalism (e.g. 
Efstathiou & Rees 1988; Haiman & Loeb 1998). The AGN 
luminosity is supposed close to the Eddington value dur- 
ing a "duty cycle" that can be as short as 10 6 to some 10 7 
years (e.g. Haiman & Loeb 1998; Martini & Weinberg 2001; 
Wyithe & Loeb 2002) or longer, of the order of some 10 s 



years (e.g. Haehnelt & Rees 1993; Haehnelt et al. 1998). 
With time, these models have become more and more com- 
plex and complete. The accretion rate onto the MBH varies 
depending on the model: it can be constant or redshift- 
dependent (KaufFmann & Haehnelt 2000; Haiman & Menou 
2000). The AGN activity in many cases appears to be a re- 
current event (Siemiginowska & Elvis 1997; Hatziminaoglou 
et al. 2001). The diverging opinions as to this recurrence 
have converged, accounting for tidal interactions and merg- 
ing effects. In addition, the Eddington-limited capability of 
MBHs to accrete mass and the progressive consumption of 
the hot gas available in haloes are also expected to play an 
important role in the evolution of AGN luminosities (Cava- 
liere & Vittorini 2000). 

In several recent papers (Silk & Rees 1998; Monaco et 
al. 2000; Granato et al. 2001; Archibald et al. 2002; Volon- 
teri, Haardt & Madau 2003 and references therein), the joint 
modelling of star formation and the growth of MBHs lying in 
galactic nuclei has been considered in order to account for 
the observed correlation between the masses of the MBH 
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and the hosting spheroid (Wandel 1999; Woo & Urry 2002). 
Spheroids grow through galaxy mergers, giving rise to new 
ellipticals, and through the mass transfer from unstable discs 
(possibly not even settled down) to bulges. Although none of 
these processes is directly associated with the Major Merger 
(MM) of dark-matter haloes, some concatenation of events 
is believed to correlate both phenomena. Firstly, MMs mark 
the beginning of the orbital decay, through dynamical fric- 
tion, of the more massive galaxies of the progenitors, which 
eventually merge at the center of the new halo. Secondly, the 
shock-heated gas is redistributed and begins to cool from 
the central higher density region at the strongest rate ever 
reached during the whole halo life. Large amounts of gas 
then fall into the central galaxy making its disc easily be- 
come unstable (particularly at high redshifts where discs 
are more compact) and transfer mass into the bulge. As the 
spheroid grows some mass may reach the nuclear region and 
feed the central MBH. We therefore expect the growth of 
spheroids and MBHs not only go in parallel, but also be 
punctuated by the rate at which dark-matter haloes form 
as a consequence of MMs (Kauffmann & Haehnelt 2000; 
Haiman & Menou 2000). 

It is difficult to tell at which extent this scenario is able 
to reproduce the observed properties of AGN and their de- 
pendence with redshift. The models found in the literature 
only address partial aspects of it and often include poorly 
justified assumptions and rough estimates of the MBH host 
density. Nowadays, the development of a complete model of 
this kind has become feasible. In particular, accurate ana- 
lytical expressions for the halo MM rate are available (e.g. 
Percival, Miller & Peacock 2000; Raig, Gonzalez-Casado & 
Salvador-Sole 2001) that can be easily implemented. The 
aim of the present paper is to build such a complete, physi- 
cally motivated model and check the validity of the previous 
scenario by comparing the predicted luminosity function of 
quasars (QLF) with the most recent publicly available ob- 
servations. 

Our model is based on the three following ingredients: 
1) an accurate expression, in terms of halo mass and cosmic 
time, for the halo MM rate (Sec. 2) derived in the frame of 
the extended Press-Schechter model; 2) the evolving prop- 
erties of MBHs (Sec. 3), obtained by solving analytically 
the governing set of differential equations established un- 
der a few simple assumptions consistent with the parallel 
growth of MBHs and their hosting galaxies from the hot gas 
available in the halo; and 3) the typical lightcurves of AGN 
associated with MBHs of different masses at different cos- 
mic times, inferred under the assumption of an Eddington- 
regulated accretion rate (Sec. 4). By comparing the resulting 
theoretical QLFs (Sec. 5) with the ones observed at differ- 
ent redshifts, the best values of the free parameters entering 
the model are inferred. This leads to a fully consistent an- 
alytical expression of the QLF (Sec. 6) that recovers both 
the observed QLF shape and its redshift evolution at in- 
termediate and high redshifts. In the last section, we dis- 
cuss the possible shortcomings of the model and the impli- 
cations of our results on the process of galaxy formation 
and evolution. Throughout the present work, we assume a 
ACDM cosmology, defined by the values of the Hubble con- 
stant (in km s _1 Mpc -1 ), density parameter, cosmological 
constant, baryonic fraction, and normalisation of the power 



spectrum respectively given by (Ho, ^m, ^a, Qbh ,as) = 
(70, 0.3, 0.7, 0.02, 0.75), with h = Ho/100. 



2 HALO MAJOR MERGER RATE 

The Modified Press-Schechter (MPS) model developed 
by Salvador-Sole, Solanes & Manrique (1998) and Raig, 
Gonzalez-Casado & Salvador-Sole (1998, 2001) is a variant 
of the extended PS model (Bower 1991; Bond et al. 1991; 
Lacey & Cole 1993) that provides an unambiguous defini- 
tion of the halo formation and destruction times and other 
quantities, such as the halo MM rate, not available in the 
original extended PS model. 

By definition, haloes are destroyed in mergers yielding 
a fractional mass increase above some fixed threshold A m ; 
otherwise they survive. The formation of new haloes is then 
consistently defined by those mergers in which all partners 
are destroyed and the whole system rearranges. These are 
the merger events we refer here to as MMs. As shown by 
Raig et al. (2001), MMs are essentially binary and involve 
similarly massive haloes. More specifically, the halo inter- 
nal structure found in high-resolution JV-body simulations 
(Navarro, Frenk & White 1997) is recovered (Manrique et 
al. 2003) provided A m ~ 0.5. This is the value of A m we use 
hereafter. 

In this model, the comoving number density of haloes 
with masses in the range [Mh.Mh + cIMh] at a cosmic time 
t is given by the usual PS mass function 



/(M H ,t) dM H = 



2 po 6 c (t) 



ttM h ct(Mh) 



dln<T(M H ) 



dlnM H 



exp 



s 2 c (t) 



2a 2 (M H ) 
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where S c (t) is the linear extrapolation to the present time 
to of the critical over-density for collapse at t, a(Mn) de- 
notes the r.m.s. mass fluctuation of the density field at to 
smoothed over spheres of mass Mh, and po is the current 
mean density of the universe. On the other hand, the merger 
rate at t of progenitors with mass Mp giving rise to haloes 
with masses in the range [Mh, Mh +(1Mh] is (Lacey & Cole 
1993) 
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Note that the previous expression corresponds to a transi- 
tion rate, i.e. it is defined per halo of mass Mh. To obtain 
the increase per unit comoving volume of the total number 
of haloes with masses [Mh, Mh + cIMh] owing to mergers 
among progenitors with masses [Mp,Mp + dMp] one must 
take R(M P — > Mh , i) dM H times f(M P ,t) dM P . 

The required total MM rate per unit comoving volume 
giving rise to new haloes with masses [Mh, Mh + cIMh] is 
therefore the integral of the latter quantity over the mass of 
the primary progenitor in the appropriate interval 

fl M M(M H ,t)dM H = 
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Figure 1. MM rate integrated over one 0.5-dex interval in halo mass (left panel) and typical interval between two MMs (right panel) 
vs. redshift, for three halo masses: 3 X 10 11 (dashed line), 3 X 10 12 (dotted line) and 3 X 10 13 M Q (solid line). For the integrated MM 
rate, these masses correspond to the lower bounds of the respective integration intervals. 



dM H 



M H /(1+A m ) 



R(M P -> Mn,t)f(Mp,t)dM P . (3) 



M H /2 



The behaviour of such a MM rate is shown in the left panel of 
Fig. 1. The associated formation rate, -Rform, defined again 
as a transition rate, is then equal to the MM rate given 
by equation (3) divided by f{Mn,i)AMn. The goodness of 
the preceding expressions has been checked against iV-body 
simulations by Raig et al. (2001). 

We want to outline the difference between the AGN 
host abundance resulting from the halo MM rate derived 
here, equation (3), and that considered in the vast majority 
of preceding models, simply given by the PS mass function, 
equation (1). The PS mass function implies a positive evolu- 
tion (i.e. it increases with increasing redshift) of the comov- 
ing number density of haloes with masses below the typical 
mass for collapse. However, this positive evolution is not 
enough for the predicted QLF to reproduce the strong pos- 
itive evolution of the observed QLF. As shown in Appendix 
A, the halo MM rate essentially behaves as the product of 
two PS mass functions. This boosts the positive evolution 
of the predicted QLF, which is decisive in making it match 
the observations. 



3 EVOLVING PROPERTIES OF MBHS 

The AGN luminosity depends on the mass Mbh of the as- 
sociated MBHs and the amount of mass AMbh they ac- 
crete (see Sec. 4). In our model, the secular evolution of 
the mass of a MBH corresponds to a sequence of short 
episodes of intense accretion, with characteristic timescale 
r acc each, between consecutive MMs typically separated by 
tmm = 1 /-Rform- Since r acc , essentially the duty cycle of 



AGN, is expected to be of the order of 10 6 — 10 9 yr, we 
have Tacc <C tmm (for a redshift z smaller than ~ 6, see the 
right panel of Fig. 1). Except for the T aC c interval after each 
MM, the accretion rate is thus negligible and the AGN is 
dormant. 

Given the random character of MMs, haloes with the 
same initial mass will evolve in different ways. Similarly, 
MBHs located in haloes of the same mass Mh will be found 
at different growth phases. Here we are not interested in 
modelling the stair-like (Lagrangian) growth of any individ- 
ual halo or any individual MBH, but the smoother (Eule- 
rian) evolution of the average mass at t, Mh(£), of haloes 
with the same initial mass and of the average properties at t 
of MBHs, namely M B H(M H ,t) and AM BH (M H ,() in haloes 
of mass Mh- 

According to the definition of tmm , the Lagrangian evo- 
lution of the masses Mh and Mbh of objects surviving at t 
approximately satisfy 



Mh (t) 



dM H _ 

dt ~ T U M[M H (t),t] 

dM BH _ AMbh [M H 



dt 



tmm [M H (*),*] 



(4) 
(5) 



As shown in Appendix B, one can solve analytically equa- 
tion (4) to obtain the evolution of Mh under the effects 
of MMs (we are neglecting the growth of haloes through 
accretion). To solve equation (5), however, we need to 
relate AM B h with M B h- According to Cavaliere & Vit- 
torini (2000), we must distinguish between the Eddington- 
limited or the supply-limited accretion regimes, depending 
on whether the amount of material afforded by the MBH 
radiating close to the Eddington limit (see Sec. 4.1) is larger 
or smaller than the one that can be accreted, respectively. 
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Figure 2. Lagrangian (left panel) and Eulcrian (right panel) evolution of A/bhj the instantaneous average mass (dotted line), and 
AMbhi the average mass increase per MM (dashed line), of MBHs located within dark- matter haloes with evolving or fixed mass Mg 
(solid line) equal to 10 13 M© at z c = 2 where the change from the Eddigton-limitcd regime to the supply-limited regime takes place. 
The values of the remaining parameters k, eq and e g are those quoted in the third column of Table 1. 



As the Eddington luminosity is proportional to the mass 
of the MBH (see eq. [14]) and the duty cycle can be assumed, 
in a first approximation, with a fixed typical value, in the 
Eddington- limited accretion regime, AM B H(M H ,i) should 
be essentially proportional to Mbh(Mh,£), 



AM BH (M H ,() = k Mbh (Mh, i) . 



(6) 



As shown in Appendix B, for a fixed value of k taken as 
a free parameter of the model and some initial conditions 
set by the matching solution in the supply-limited regime, 
the assumption (6) leads to analytical expressions for the 
wanted dependences on Mh and t of Mbh and AMbh. 

The solution corresponding to the Eddington-limited 
accretion regime is expected to hold at high redshifts where 
the amount of hot gas in haloes is large and its radiative cool- 
ing very effective owing to the high inner density of haloes. 
At small enough redshifts, however, a substantial fraction of 
the hot gas in the initial haloes will be already consumed 
and the accretion of MBHs will become supply-limited. We 
then expect the mass accreted by MBHs after the MM of 
their hosting haloes to depend on the small amount of mass 
accreted by the galaxy at the same occasion. The evolution 
of the two quantities Mbh(Mh, t) and AMbh(Mh, t) is then 
inseparable from that of the masses of the hot gas available 
in the halo, M g (Mn,t), and the central galaxy, Mg(Mh,<). 

What do we know about the relation between Mbh 
and Mg? As suggested by Kormendy & Richstone (1995), 
the masses of MBHs and spheroids are closely related. 
Magorrian et al. (1998) found a simple proportionality be- 
tween the MBH mass and the luminosity of elliptical galax- 
ies, while subsequent works have shown that the relation 
A^bh oc Uv' 5±0 ' 5 , with ov the bulge velocity dispersion of 



the host galaxy, gives a somewhat tighter correlation (Fer- 
rarese & Merritt 2000; Gebhardt et al. 2000; Merritt & Fer- 
rarese 2001; Tremaine et al. 2002). Combining the latter 
with the relation observed between a v and circular veloc- 
ity (Ferrarese 2002), likely more closely related to the halo 
virial mass, results in a proportionality between the MBH 
and halo masses, Mbh and Mh, which is also supported the- 
oretically (see Ferrarese 2002 and references therein). This 
latter proportionality might be used to infer the evolution 
with cosmic time of Mbh, during the supply- limited accre- 
tion regime, from that of Mh valid at both the Eddington- 
limited and supply-limited accretion regimes. All the previ- 
ous observations refer, however, to nearby galaxies and any 
extrapolation to higher redshifts would be very risky. On the 
other hand, a definite time-evolution of both Mbh/Mh and 
Mg /Mh automatically follows, in the supply-limited accre- 
tion regime, from the assumed parallel growth of Mbh and 
Mg expected to occur as a consequence of halo MMs. 

This should take place according to the following simple 
assumptions or approximations, i) The baryonic component 
of a halo consists essentially of hot gas and one main central 
galaxy. Note that, if the two central galaxies in the progeni- 
tor haloes merge, the new central galaxy will automatically 
include their baryonic masses. In contrast, we neglect those 
baryons locked in surviving low mass satellites, ii) The av- 
erage mass accreted by the central galaxy after a MM is 
roughly a constant fraction of the hot gas mass M g avail- 
able at the moment of the MM. iii) The average mass AMbh 
accreted by the MBH after some delay that can be also ne- 
glected in front of tmm is a constant fraction of the mass 
accreted by the hosting galaxy as a consequence of the MM 
or, equivalently, due to assumption (ii), a constant fraction 
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e g of M g at the time of the MM triggering the whole process. 
Besides, assumptions (ii) and (iii) lead to a constant ratio 
eg between the average mass Mbh of the MBH and the av- 
erage mass Mg of the central galaxy, as roughly observed. 
We can therefore write 



M g (M H , t ) + M G (Mh , t ) = M b (M H ) 
Mbh (Mh, t) = e G M G (M H ,i) 
AM BH (M H ,() = e g M s (M H ,t) , 



(7) 
(8) 
(9) 



where Mb is the average baryonic mass in haloes of mass 
Mh, equal to the constant ratio /b = Qb/^m times Mh, 
with Qb and n m the baryonic and total density parameters, 
respectively, and eg and e g are two free parameters of the 
model. 

The dependences on Mh and t of all previous aver- 
age quantities, in particular of Mbh and AMbh, in the 
supply-limited accretion regime are derived in Appendix B 
by solving analytically the set of differential equations im- 
plied by the assumptions (7)-(9) and equations (4) and (5). 
As an initial condition it is assumed that, at some redshift 
z c marking the change between the Eddington-limited and 
the supply-limited accretion regimes, all haloes have essen- 
tially the same baryonic fraction in the form of hot gas, 
7 = Mg/Mb = 7c- The redshift z c is taken as another free 
parameter of the model. In contrast, the match of Mbh and 
AMbh in the two regimes at z c leads to a one-to-one corre- 
spondence between 7 C and k and there is no need to intro- 
duce any extra freedom in the model (see Appendix B). 

The predicted Lagrangian evolution of Mh, Mbh and 
AMbh is shown in the left panel of Fig. 2. In the plot, we 
represent the results obtained for haloes of 10 13 M at z c , 
although all halo masses lead to very similar results. As can 
be seen, AMbh and Mbh, which initially differ by the con- 
stant factor k, become equal shortly after z c . In the latter 
supply-limited accretion regime, Mbh is also equal to eG Mg 
(by assumption [8] ) and approximates to eg M b as t tends to 
infinity (where all baryons tend to be locked in the central 
galaxy). On the other hand, the trend with time of the ratios 
Mbh/Mh and AMbh/Mh changes from strongly increasing 
to slightly decreasing in the passage from the Eddington- 
limited to the supply-limited accretion regimes as a conse- 
quence of the coupling, in the latter, of AMbh with the de- 
creasing gas mass fraction. But what determines the evolu- 
tion of the AGN luminosity is rather the Eulerian evolution 
of Mbh and AMbh (see Sec. 4), which is represented in the 
right panel of Fig. 2. This plot shows that both quantities 
increase with increasing time until z c , where the values of 
Mbh begin to slightly decrease converging swiftly to those 
of AMbh, forming a plateau for z below ~ 1.4. 



4 AGN TYPICAL LIGHTCURVES 
4.1 Accretion episode 

Let us now focus on the typical evolution, during one in- 
dividual accretion episode, of a MBH with initial mass, at 
time t, equal to the average mass of MBHs in haloes of mass 
Mh at that moment, accreting a total mass of gas also equal 
to the average value. In order to clearly distinguish between 
the initial mass of the MBH, given by MBn(Mn,t) (vary- 
ing in a timescale tmm), and the mass of the MBH during 



the accretion episode (varying in a characteristic timescale 
Tacc <S tmm), we denote this latter simply by M, without 
subindex BH. 

As mentioned, some time after a MM of haloes we ex- 
pect radial infall of material to start feeding the MBH at the 
centre of the new halo. This central MBH can be the result 
of the merger of the two MBHs in the respective progenitor 
haloes, since the respective central galaxies can merge and 
their central MBHs may eventually fall into the centre of 
the new spheroid (also due to dynamical friction within the 
new galaxy) and merge. But only one of these two MBHs 
may finally be found at the nucleus of the central galaxy. 
In any event, we do not need to consider the frequency of 
MBH mergers since the initial mass assumed for the accret- 
ing MBH, equal to the average mass of MBHs in haloes with 
Mh at t, implicitly accounts for it. 

Material is expected to be supplied to the central MBH, 
initially in an increasing phase, then in a decreasing one, 
during a time of the order of r acc . Therefore, we assume 
that the accretion rate M follows a bell-shaped form 

m'{0) 



M(r) = AM BH (M H ,i) 



(10) 



determined by the dimensionless phenomenological function 



m'{6) = - exp (-0) [1 - exp (-361)] 



(11) 



with t the time elapsed since the beginning of accretion t, 
6 = r/r acc , and a dot denoting differentiation over r. The 
typical evolution of the mass M of the MBH during such an 
accretion episode then follows by integration of equations 
(10) and (11) 

M(r ) = Mbh(Mh, t) + AMbh(Mh, t) m{0) , (12) 
with m(9) equal to the integral of m! from to 9 

m(0) = 1 - | exp(-0) + | cxp(-40) . (13) 

AMbh and Mbh in equations (10) and (12) are the 
evolving average properties of MBHs introduced in Section 
3. The value of T acc , which should be related to the dynam- 
ical properties of the MBH close environment, is difficult to 
tell. For this reason, we take r acc as another free parameter 
of the model. 

4.2 Radiation model 

The accretion of gas onto the MBH reactivates the AGN. 
An accurate model of the resulting lightcurve would re- 
quire a relativistic approach. But this is outside the scope of 
the present study, aimed merely at understanding the main 
trends of the observed QLF and its redshift evolution. For 
this reason we shall neglect any relativistic correction. 

Under the classical black hole paradigm, a source with 
mass M powered by spherical accretion can reach a maxi- 
mum luminosity, called Eddington luminosity, given by 

47rcGpaccM 



Lsdd = 



(14) 



OThCEddWe 

with G the gravitational constant, c the speed of light, p acc 
and n e the mass density and electron density, respectively, 
of the accreted matter, <r T h the Thompson cross-section and 
C*Edd — 0.2 a dimensionless constant which accounts for 
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Figure 3. Predicted AGN lightcurves (left panel) and corresponding Eddington efficiencies (right panel) during an accretion episode 
onto the average MBH within haloes of typical mass 10 13 M at redshifts equal to 3.0 (full line), 2.6 (dotted line), 2.2 (dashed line), 1.7 
(dot-dashed line) and 0.3 (triple dot-dashed line). In the case of £Edd we on ly plot the solutions for z ^ z c since the remaining curves 
overlap with that of z c = 2. The values of all the parameters used are given in the third column of Table 1. 



the variations with electron energy of the Compton cross- 
section, averaged with a standard quasar SED (Ciotti & 
Ostriker 2001). 

If, at some radius r at a given time r, the luminosity 
exceeds the Eddington limit, the radiative acceleration will 
exceed the gravitational acceleration and accretion will stop. 
This will cause the radiation to stop soon after, allowing 
accretion to start again, and so on. In this way, accretion 
and radiation will oscillate on a characteristic timescale of 
the order of the infall time at r. In fact, provided that this 
characteristic time is short enough, the luminosity of the 
AGN will tend to adjust to the accretion rate set by external 
conditions linked to the MM event, thus reaching a (possibly 
time-averaged) quasi-stationary regime. The characteristic 
timescales of light variations observed in AGN (years or less) 
are considerably shorter than the timescale of accretion onto 
the MBH (at least 10 6 yr). Hence, such a quasi-stationary 
regime should be reached swiftly during accretion. 

In these circumstances, the radiative pressure at any 
radius r results in a reduced effective gravity 



ffcft(r) = 



EEdd ) 



(15) 



which we have expressed in terms of the Eddington efficiency 
tEdd = i/I/Edd- The bulk of the energy release comes from 
the region with maximum accretion efficiency, which occurs 
at a radius 



^last — &last 



2 GM 



(16) 



corresponding to the last marginally stable Keplerian orbit 
(Chakrabarti 2000), with 2GM/c 2 the Schwarzschild radius 
and fci aB t a proportionality factor typically ranging between 



one half and two or three, depending on the form of the 
orbit and the angular momentum (hence, the metric, either 
Schwarzschild or Kerr) and the topology (e.g. simple, binary, 
toroidal) of the MBH. In the present calculations, we do not 
distinguish however among all these possible configurations 
and adopt the intermediate value of feiast = 1-5. According to 
equation (15), the terminal infall velocity «i ast for accreted 
material is, therefore, 



2 

Wlast 



2^(1 

Hast 



^Edd) ■ 



(17) 



Since the bolometric luminosity is the kinetic energy of the 
accreted material converted into radiation, we have 

L=^MvLt, (18) 

and, by substituting equations (16) and (17) into equation 
(18) and dividing it by L^dd, we find 



£Edd = 



1 + 



S-K fclgst Gp accM 

c u T h C*Edd n e M 



(19) 



Finally, by taking into account the definition of the Edding- 
ton efficiency tEdd given above, we are led to the following 
non-linear expression for the AGN bolometric luminosity 



hi 

L 



M 8 M s 

M m 



(20) 



in terms of the two functions, M and M, defining the accre- 
tion episode (eqs. [10] - [13]) and the three constants 
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The resulting typical AGN lightcurves, L(t), and the 
corresponding Eddington efficiencies, e E dd(i~) = L(t)/Ls x 
Ms/M(t), are shown in Fig. 3 for the same fixed halo mass 
plotted in Fig. 2 and several illustrative redshifts. Note that, 
in our model, eEdd is not constant along the duty cycle but, 
like the bolometric luminosity, it is a bell-shaped function 
of time taking values substantially lower than unity. The 
maxima of both L and eEdd depend on the Eulerian evo- 
lution of Mbh and AMbh (see the right panel of Fig. 2) 
through equations (10) and (12). From equation (19), we 
have that eEdd is constant as long as AMbh is proportional 
to Mbh , namely, during the whole Eddington regime and at 
redshifts below ~ 1.4. Between these two phases, right after 
z c , it diminishes with increasing time. On the other hand, 
from equation (20) we have that, as long as AMbh is pro- 
portional to Mbh, L is also proportional to them. There- 
fore, the luminosity increases with time during the whole 
Eddington- limited accretion regime and diminishes, after z c , 
to a constant value where Mbh and Abh become equal. This 
implies a negative luminosity evolution of AGN associated 
with haloes of a fixed mass during the Eddington-limited 
accretion regime, a positive evolution in the initial phase of 
the supply- limited accretion regime until z ~ 1.4 and a null 
evolution at lower redshifts. It is worth noting that the lu- 
minosity also increases with increasing halo mass, but not 
the Eddington efficiency. 



5 AGN LUMINOSITY FUNCTION 

The (differential) AGN luminosity function, <p(L,t), is de- 
fined through the relation 



d 2 N(L, t) = dn(L, t) dV = <j>(L, t) dL dV , 



(21) 



where N(L, t) and n(L, t) are the number and comoving 
number density, respectively, of AGN with absolute bolo- 
metric luminosity L at a given time t and dV is the element 
of comoving volume. Had all AGN the same typical strictly 
increasing (or decreasing) lightcurve L(r), we would have 



4>(L,t) = 



dn(L,t) 



dL 



± lim — — 

al->o AL 



± lim — 

A£-»0 AL 



[ R IT (t') dt' - [ flrr(t') dt' 
Jt-T(L) Jt-rlL+AL) 



± lim R TI (t) 

A£->0 



t(L+AL) 

R lr (t')dt' p 

t(L) 

t(L + AL) — t(L) 



AL 



= ±Rrr(t) 



dr 
dL' 



(22) 



where R IT is the AGN reactivation rate, the + ( — ) sign is for 
a strictly increasing (decreasing) lightcurve and the two inte- 
grals on the right of the first equality give the comoving num- 
ber density of AGN with luminosity larger (smaller) than L 
and L + AL, respectively, at t. Note that, in the derivation of 
equation (22), we have taken into account the fact that the 
characteristic accretion time T aC c is much shorter than the 
characteristic time between MMs tmm, so that, along the 
time interval of integration, R rl is approximately constant 
and equal to R TI (t). 



Therefore, in the more realistic case of AGN with typi- 
cal lightcurves dependent on halo mass Mh, we can write 



4>(L,t) 



f 

Jo 



dMn -Rrr (Mh , t) 



dn 
dL 



dT2 

dL 



(L,M H ,t), (23) 



where R IT (Mn, t) dMn is the reactivation rate of AGN in 
haloes with masses in the range [Mh, Mn+dMn] and where 
we have accounted for the fact that, as shown in Fig. 3, each 
typical quasar lightcurve has two branches, one increasing 
and the other decreasing with increasing r (indexes 1 and 
2, respectively). 

Hence, to infer the QLF predicted by our model at any 
time or redshift, we must simply write, in equation (23), the 
quasar reactivation rate in terms of the halo MM rate, 



R IT (M H ,t) = v -Rmm(Mh , t) . 



(24) 



with Emm given in equation (3), and take the functions 
n(L, Mu,t) as the inverse of the increasing and decreas- 
ing parts of the typical luminosity curve L(t,Mh,<) given 
in equation (20). The proportionality factor v in equation 
(24) accounts for possible inaccuracies in the assumptions 
adopted as well as for possible biases affecting the observa- 
tional data. For example, although our model presumes one 
single quasar event per halo MM, this might deviate from 
reality (factor A^qso/^mm in equation [25]) as there are in- 
dications that haloes may harbour more than one quasar 
(Mortlock et al. 1999; see also Croom et al. 2001b). On 
the other hand, there are also reasons to introduce such 
a normalisation factor from the observational point of view. 
More specifically, no sample is ever free of incompleteness ef- 
fects and selection biases (factors C and -/VQso,obs/-'VQso,tot 
in equation [25], respectively). Last but not least, one must 
take into account, in addition, the orientation of detected 
quasars (factor 2f2 e m/47r in equation [25]). All these ef- 
fects should not present a strong evolution in the range 
< z < 2.5, but be rather constant (with eventual variations 
in some specific redshift ranges due to redshift-dependent 
observational biases). As several of these effects are rather 
uncertain, we take 



-C x 



2Q, 



em ■"QSO.obs Aqso 



4-7T A^Qso.tot Nmm 
as another free parameter of the model. 



(25) 



6 COMPARISON TO OBSERVATIONS 

As explained in Appendix B, we do not expect z c to sub- 
stantially vary with halo mass, but the situation is less clear 
for the other quantities taken as constant parameters in the 
present model. For instance, by combining the scaling rela- 
tion (A4) with the observed Mbh — cr v correlation, we could 
infer an expression for eG as a function of halo mass and 
redshift. However, letting these quantities to depend on Mh 
and t or on any of these variables alone would greatly com- 
plicate the task of obtaining an analytical model and would 
introduce too much freedom in the problem. Thus, we prefer 
to assume them having fixed values in a first approximation. 
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Figure 4. Predicted QLFs (lines) compared to the observed blue QLFs drawn from the 2QZ sample (symbols) at various rcdshifts. The 
tiny error bars affecting the observational data are only statistical. Left panel: rcdshifts in the [0.3, 2.3] range; quasars are binned around 
(z}=0.435, 0.730, 1.080, 1.808, and 2.020 (bottom to top). Right panel: restriction to the two bins with (2) larger than unity where the 
model gives the best fit. The values of the parameters used are given in Table 1, columns 2 and 3, respectively. For clarity, the QLFs are 
shifted upwards by 0.2 units with each increasing (z). 



6.1 Constraints to the parameters 

So far we have implicitly assumed that all haloes, irrespec- 
tive of their mass, harbour MBHs evolving according to the 
model given in Section 3 and that all MBHs radiate accord- 
ing to the model described in Section 4. But this may not 
be the case. Only haloes in some finite mass range are likely 
to feed MBHs at the suited rate and reach the appropriate 
radiative efficiency (i.e. the capability of converting the ki- 
netic energy of the matter falling at ri ast into radiation) . We 
should therefore also examine the effects of adopting a lim- 
ited halo mass range, [M m i n , M max ]. In summary, the free 
parameters of the model are: z c , eg, £g, 7~acc, k, M min , M max 
and v. 

To find the best values of these parameters we have fol- 
lowed an iterative ^-minimisation fitting procedure to the 
observed QLF, with each free parameter of the model tak- 
ing values from a pre-constructed grid within certain initial 
intervals (see below). For the following iterations, however, 
the width of these intervals is progressively reduced and their 
centres are re-defined using the best-fitting parameter values 
of the former step in such a way that the free parameters 
are allowed to take values outside their initial intervals. (Ac- 
tually, the only parameters fitted in this way are z c , eg, 
T acc , k, Mmi n and M max ; the value of v is drawn, at each 
step, by taking the best simultaneous normalisation of the 
predicted QLFs.) 

The redshift z c marking the change from the Eddington- 
limited to the supply-limited accretion regimes is initially 
allowed to take values in the range [0,10]. The ratio be- 
tween the masses of MBHs and their central galaxies, ec, 
is allowed to vary in the range [5.5 x 10~ 4 ,6.7 x 10~ 3 ] 



indicated by observations (Kormendy & Gebhardt 2001; 
McLure & Dunlop 2001; Wandel 2002; Falomo et al. 2002; 
Wu et al. 2002). Contrarily to eg, the fraction of gas in 
the host galaxy that is accreted by the MBH at the oc- 
casion of a MM, e g , is poorly constrained. So we leave e g 
completely free. The timescale r acc is supposed to be in the 
range [10 6 ,10 9 ] yr bracketing the values mentioned in the 
literature for the duty cycle of AGN. The previous bounds 
for r acc together with equation (19) under the approxima- 
tion M = AM/(2r acc ) = fcM/(2r acc ) lead to a fractional 
mass increase, k, of the MBH during an Eddington-limited 
accretion episode in the range [0.05,50]. M m i n and M max 
are assumed to lie initially in the ranges [10 11 ,10 1 ] M© 
and [10 13 , 10 14 ' 5 ] Mq, respectively. Finally, the normalisa- 
tion factor v is also left completely free. 



6.2 Results for the 2dF QSO sample 

Quasars are not always easy to distinguish from other types 
of objects. Moreover, the different samples gathered in differ- 
ent wavelengths are usually small. In this concern, the 2dF 
quasar redshift (2QZ) sample (Croom et al. 2001a) seems 
the most appropriate for our purposes. The 2QZ sample is 
very large, uniformly selected from well-defined criteria and 
of high overall completeness C ~ 93% to the limiting mag- 
nitude b.j = 20.8, after corrections for incompleteness (see 
http://www.2dfquasar.org/). The subsample of 2QZ used 
here comprises some 10000 colour-selected, moderately faint 
(Mb < 21) quasars, covering the redshift range [0.3,2.3]. The 
candidates were all point-like objects selected by their UV 
excess on a (U — B) x (B — R) plane, introducing an abrupt 
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Figure 5. Predicted and observed QLFs for quasars in the very 
high-z range [3,4.7] binned around (z) = 3.75 (solid line and 
squares), 4.15 (dotted line and circles), 4.3 (dashed line and di- 
amonds), 4.4 (dot-dashed lines and crosses) and 4.7 (triple dot- 
dashed lines and triangles). The values of the parameters used 
are given in column 4 of Table 1. Both predictions and data are 
shifted downwards by 0.3 units with each increasing (z). Note 
that the two first points on the left are only upper limits. 



decrease of the efficiency of the selection both at z > 2.3 
(UVX criterion) and at z < 0.3 (due to morphological clas- 
sification). To perform the fit, the corresponding QLFs for 
several redshift bins derived by Croom et al. (2001a) assum- 
ing Ho = 50 km s _1 Mpc -1 , Q, m — 1 and Qa = have been 
converted (in luminosities and densities) to the ACDM cos- 
mology used here and the bolometric luminosities predicted 
by our model transformed to absolute luminosities in the 
blue band by considering a constant ratio Lb/L Bo i = 0.023, 
inferred from a typical quasar SED (Elvis et al. 1994). 

In the left panel of Fig. 4, we plot the predicted and 
observed QLFs for the best-fitting values of the parameters, 
listed in the second column of Table 1 . In the third column of 
the Table, we also quote the best values of the parameters for 
the data restricted to 1.27 < z < 2.3, where the model gives 
a much better fit (see the right panel of Fig. 4) . The possible 
origin of the poor fit at low redshifts is discussed in Section 
7. As can be seen from Table 1, all parameters take values 
in the respective expected ranges. The value of M m i n quoted 
in the second column is, however, essentially an upper limit 
since the lower error bound is not constrained by the data. 
This is because varying this lower mass boundary essentially 
modifies the fainter part of the luminosity function while, 
at all redshifts, the cut-off due to the limiting magnitude 
bj = 20.8 is effective at magnitudes brighter than the cut- 
off due to M m i n . Besides, at low z there is another cut-off at 
Mb = —23.0 due to the selection of stellar morphologies. A 
more realistic model, likely giving better fits to the observed 
QLFs at the faint end, would require a progressive decline 



of the integrand in equation (23) beyond some halo mass 
instead of a sharp cut-off at M m i n , although at the cost of 
increasing the number of free parameters. In this concern, 
the best value of M m i n should be taken with caution. 

The factor A?Qso,obs/A'QSO,tot appearing in the decom- 
position of v given in equation (25) reflects two combined 
effects: the intrinsic properties of the objects and selection 
biases. The 2QZ quasars are blue objects with redshifts lower 
than ^2.3, selected by their UV excess. Therefore, this fac- 
tor represents here the number of low to intermediate blue 
quasars over the total quasar population. Several recent re- 
sults based on X-ray observations (e.g. Rosati et al. 2002) 
show that the reddened, or type-II, quasars could be as many 
as the blue ones detected in the optical passbands and, there- 
fore, UVX quasars could be only a fraction (but a large one) 
of the overall optical quasar population, estimated to be of 
the order of 50—90%. The factor 2Q cu1 /4tv, accounting for 
the orientation effects, is ~ 0.25 according to the unified 
models (see e.g. Elvis 2002). As mentioned in Section 5, 
JVqso/^Vmm ~ 1, although (slightly) greater values are the- 
oretically possible. Accordingly, a most crude estimate of 
the value of v gives an upper limit of 10% to within a factor 
of around 2. Thus, the predicted value of ~ 10% also lies 
within the expected range. 



6.3 Results for QSOs at very high redshifts 

We have also compared the predictions of our model to the 
very high redshift QLF, compiled from the points of Fan et 
al. (2001), Kennefick et al. (1995), Kennefick et al. (1996) 
and Schmidt et al. (1995). The best fit is once again sat- 
isfactory (see Fig. 5), although the values of some param- 
eters move substantially apart from those previously found 
(compare column 4 with column 3 of Table 1). However, 
taking into account the large uncertainty affecting the new 
values, we see that such differences are not significant. In 
fact, the results of this high-z fit must be taken with cau- 
tion because of possible systematic effects, such as lensing, 
affecting the data (see e.g. Wyithe & Loeb 2002). Moreover, 
the small error bars associated with the values of the param- 
eters drawn from the 2QZ sample are likely underestimated 
since they only account, once again, for statistical errors in 
the observed QLFs, almost negligible in this case owing to 
the large number of quasars included in each bin. The only 
parameters that, despite all, seem to change are r acc and, 
to a lesser extent, M m i n . As mentioned, we should not at- 
tach too much importance to the exact value of this latter 
parameter as a smooth cut-off in halo mass is likely more 
realistic than a sharp one. Concerning T acc , the tendency it 
shows to decrease with increasing z was foreseeable since all 
physical sizes and separations become smaller and, hence, 
all dynamical processes go faster at higher redshifts. In par- 
ticular, the typical interval between MMs diminishes with 
increasing z (see Fig. 1), thus r acc must also do so if the 
condition r acc <C tmm is to be satisfied. 



6.4 Results for the integrated quasar light 

An alternate way of estimating the best values of the param- 
eters is matching the quasar light density over the redshift 
range [0,6]. Note, however, that the model is very poorly 
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Figure 6. Integrated quasar light density vs. rcdshift compared to data in optical (Shaver et al. 1996) and X-rays (Miyaji et al. 2000), 
and the ionising flux estimated from the absorption spectra of high- 2 quasars (McDonald & Miralda-Escudc 2001). The model predictions 
are computed for the set of parameters shown in Table 1, column 5. Data and model predictions are normalised to z = 2. 



constrained by the observations at z > 4 owing to the large 
fraction of the quasar population that is likely to be dust- 
enshrouded (Xia et al. 2001; Rosati et al. 2002; Norman et 
al. 2002). Fig. 6 compares the model predictions with avail- 
able data (the gray-shaded region). The lower limit of the 
observations arises from the optical data by Shaver et al. 
(1996), while its upper limit is deduced from X-ray mea- 
surements (e.g. Miyaji et al. 2000). We also include the data 
points coming from an estimate of the photoionisation rate 
based on the observed absorption spectra in high-z quasars 
(Fig. 3 of McDonald & Miralda-Escude 2001). In fact, the 
data used to adjust our model are the latter (normalised at 
2 = 2) since they are better adapted to our \ 2 fitting pro- 
cedure. Caution must be paid, however, to the fact that the 
total ionizing flux is believed to be dominated, at z > 4, by 
that of the blue stellar population. 

The predicted behaviour of the quasar light density is 
compatible with the observed one for intermediate and high 
redhifts, while there is, once again, less evolution than ob- 
served at z < 1. Interestingly enough, the best values of 
the parameters obtained from the fit to the estimated ion- 
izing background, given in column 5 of Table 1, are in gen- 
eral close to and, in any case, compatible with those ob- 
tained from the 2QZ sample at z > 1. The only excep- 
tions are r aC c and M m i n . The former takes, not surprisingly, 
an intermediate value between those quoted in columns 3 
and 4. On the other hand, the slightly smaller value of 
the latter is even more consistent than the previous ones 
with the lower boundary of the observed mass range of cur- 
rent MBH masses, 10 6 < Mbh/Mq < 10 10 , which corre- 
sponds to 10 9 ' 6 < M H < 10 13 - 6 at the limit t -> oo where 
Mbh — > cgMg ^ eG/bAfH- This overall good agreement 
might indicate that quasars are the major responsible of the 



ionizing background at z < 4 where the fit is very good or, 
alternatively, that the ratio of the contributions to this back- 
ground of quasars and young stars remains approximately 
constant in that redshift interval. 



7 DISCUSSION 

We have developed a fully consistent analytical model of 
QLF from reasonable physical assumptions and used avail- 
able data on the QLFs at various redshifts to fix the values 
of the free parameters. The predicted QLF agrees with the 
available observations at redshifts higher than ~ 1. 

The fact that, the more accurate halo MM rate used 
in our model behaves essentially as the product of two PS 
mass functions instead of just a single one as often adopted 
explains the more marked positive evolution shown by the 
inferred QLF. Yet, the evolution so predicted is not strong 
enough to match the observations. The QLF calculated from 
the adopted MM rate but assuming a MBH mass propor- 
tional to the halo mass and the associated quasar radiating 
close to the Eddington regime fails to reproduce the ob- 
served data. A successful approach requires also modelling 
the typical evolution on cosmological timescales of MBHs 
and the typical lightcurve of AGN. In the present study, 
such a theoretical lightcurve rests on the assumptions of a 
phenomenological bell-shaped accretion rate and a physi- 
cally grounded self-regulated radiation, while the modelling 
of the MBH growth is based on reasonable assumptions em- 
anating from the expected connection between halo MMs 
and the reactivation of AGN. 

The different samples of quasars studied with redshifts 
ranging from 1.27 to 6, lead to values of the parameters 
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Table 1. Best-fitting parameters, column 1, corresponding to four data sets. In columns 2 and 3, the observed QLFs drawn from the 
2QZ survey in the two redshift intervals quoted. In column 4, the high-z QLFs drawn essentially from the SDSS survey (Fan et al. 2001) 
and, in column 5, the integrated quasar light estimated from the absorption spectra of high-z QSOs (McDonald & Miralda-Escudc 2001). 
z c is the redshift at which the accretion regime changes from Edington to supply-limited, k is the fractional mass increase of the MBH 
in the Edington-limitcd regime. €q is the mass ratio between the MBH and the host galaxy and e g the ratio between the mass accreted 
by the MBH in a MM event and the total gas mass available, both in the supply-limited accretion regime. r acc is the typical timescale of 
accretion onto the MBH (essentially the quasar duty cycle). M m j n and M max are the lower and upper masses of haloes likely to harbour 
quasars. Finally, v is the fraction of MMs giving rise to the reactivation of one observed quasar (in the case of the integrated quasar light 
there is no value for v since both the predictions and observations are normalised at z = 2). Errors correspond to la intervals, calculated 
in the standard way for a multiparamctric fit by \ 2 minimization. The lack of errors in the second column is due to the poor fit obtained 
in this case, which makes such an estimation meaningless. 
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that do not differ significantly. The only relevant exception is 
r acc , the typical timescale of accretion onto the MBH, which 
shows a clear trend to diminish with increasing z as actually 
expected. It is worth emphasizing that the values found for 
all the parameters are very reasonable. In particular, those 
of T acc , ranging from 10 to 2 x 10 yr according to z, arc 
consistent with the most recent estimates (Yu & Tremainc 
2002; Schirber & Bullock 2003). The ratio e G = M BH /M G 
is found to be ~ 2.5 x 10~ 3 , in agreement with observa- 
tions (Ferrarese & Merritt 2000; Gebhardt et al. 2000). The 
fraction e g of hot gas that ends up in the central MBH af- 
ter a MM in the supply-limited accretion regime is of a few 
percent. The value ~ 2 of the parameter k determining the 
accretion rate in the Eddington-limited regime is just what 
one would expect from equation (20) for a quasar of lumi- 
nosity Lg, representative of a MBH of 10 8 M Q , and a typical 
accretion rate of Ms during a duty cycle of 10 8 yr. On the 
other hand, according to Cavaliere & Vittorini (2000), the 
transition between the Eddington-limited and the supply- 
limited accretion regimes is expected to occur at the typi- 
cal formation of groups, at z equal to 2.5—3.0, where MMs 
of haloes of the galactic scale become rarer and the accre- 
tion onto MBHs progressively declines. In our model, such 
a transition takes place when the gas fraction in baryons 
becomes so low that the product e g M g falls below /cMbh, 
thus becoming the limiting factor. This occurs at the epoch 
where Mg/M s ~ e s /{keo), with the latter ratio taking a 
value of order unity, which is satisfied at z c ~ 2 not far from 
the redshifts preferred by Cavaliere & Vittorini (2000). The 
hosting halo masses range from ~ 10 11 M© to ~ 10 14 M©, 
implying that quasars would be found within central galax- 
ies of haloes from the galactic scale to the scale of galaxy 
groups, the structures that are currently believed to be the 
harbouring environments of the brightest quasars (Jager et 
al. 2001). Finally, the inferred fraction of observed quasars 



per halo MM is found to be v ~ 5-10%, again pretty close 
to the expectations. 

All these results, i.e. the good behaviour of the predicted 
QLFs at z > 1 and the reasonable best-fitting values of the 
parameters, give strong support to the proposed connection 
between halo MMs and AGN lightning through the triggered 
simultaneous feeding of the spheroid of the main galaxy and 
its central MBH. The fact that our model does not match 
the observed evolution of QLFs at z < 1 is also not unex- 
pected in this scheme because of the failure for late times 
of that central hypothesis. The two complementary mech- 
anisms connecting halo MMs and the reactivation of AGN 
imply the parallel growth of the host galaxy (by cooling 
flows and the merger of the main galaxies of the progenitor 
haloes). Consequently, they can be expected to operate as 
far as haloes are typically of the galactic scale, when galaxies 
are growing. When haloes reach the scale of galaxy groups, 
galaxies within them no longer grow in size, but just in num- 
ber. Therefore, none of the above mentioned processes may 
then be effective. At such late epochs, there may still be ra- 
dial flows of gas into MBHs lying at the centre of galaxies, 
but they should be triggered by other mechanisms — such as 
internal instabilities of galaxies or tidal interactions among 
galaxies within groups and clusters — uncorrelated with halo 
MMs. 

These results also indicate the necessity of dealing more 
accurately with the physics of baryons, that is, with the 
processes governing the evolution of the hot gas trapped in 
haloes, as well as with the growth of galaxies and their inter- 
actions. In other words, a QLF model spanning from high to 
small redshifts should necessarily be coupled to a detailed 
model of galaxy formation and evolution as in the recent 
work by Menci et al. (2003). There are other simplifications 
in our model that might be also improved in future works, 
although they do not seem to play a crucial role on our 
main conclusions. There is, for instance, no consideration 
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of the detailed physics at the vicinity of the MBH and our 
radiation model does not account for relativistic effects. De- 
tailed study of the physics around the MBH could provide, 
for example, a natural explanation to the decrease of the 
number of bright quasars (Nitta 1999). Likewise, the use of 
an isotropic model to describe highly anisotropic radiation 
fields around quasars is likely to introduce different kinds of 
biases. The first one, due to the requirement that the line-of- 
sight lies in the emission cone, is of purely statistical nature 
and should be reasonably accounted for by the normalisation 
of the MM rate. A second source of bias arises from the rela- 
tion between the bolometric luminosity and the MBH mass, 
through the Eddington efficiency, which assumes isotropy 
of the emission. This approximation, which cannot be cor- 
rected by a shift of the luminosity function along the mag- 
nitude axis, results in erroneous masses associated with the 
observed luminosities and, therefore, in erroneous merger 
rates. Finally, one must not forget that the comparison with 
observational data usually takes into account only statisti- 
cal errors (very small in such large samples as the 2QZ). 
Systematic effects, however, like the variations of the ability 
of redshift determination depending on the emission lines 
present inside the observable wavelenght range, are likely to 
have a major contribution (Croom et al. 2001a). 
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APPENDIX A: KINETIC APPROACH FOR 
THE HALO MM RATE 

In the simplest case of an n — 0, power-law power spectrum 
of density fluctuations (where there is no correlation among 
haloes), the merger rate is simply given, in the kinetic theory, 
by the usual expression 

R(Mn,t)dM H = 



/(M H , t) dM H f(MpH,t) dM P {v S(v))p (1 + zf , (Al) 



with Mp and Mpu = Mh — Mp the masses of the primary 
and secondary progenitor haloes, respectively, (vS(v))p the 
product of the relative velocity between the two progeni- 
tor haloes and their cross section for merger averaged over 
those velocities, and (1 + z) 3 the factor that brings the lat- 
ter product, expressed in natural units, to comoving units 
as the quantities entering in this expression. A reasonable 
approximation to adopt is 



{v S(v))p = K nrp ap . 



(A2) 



with rp and ap the gravitational radius and internal 3-D 
velocity dispersion, respectively, of the primary progenitor, 
which can be estimated from their scaling relations with Mp . 
Indeed, provided all relaxed haloes have a common density 
contrast S, one has within the radius r that M = 8r [i p = 
Sr 3 po (1 + z) 3 , with p the mean cosmic density at redshift 
z, leading to 
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Figure Al. Logarithm of the halo formation rate vs. redshift 
in the Einstein-de Sitter cosmology with an n = power-law 
spectrum. Two sets of curves have been inferred from the MPS 
model and the kinetic theory (the latter with K = 1 in equation 
[A2]), for hosting halo masses M H = 10 12 M (solid lines), M H = 
10 13 M Q , (dotted lines) and 10 14 M Q (dashed lines). On the top, 
logarithm of the corresponding ratio, which is constant and equal 
to 1550 within 2% out to z = 20 for M H = 10 12 M and within 
12% out to z = 3 for M H = 10 14 M . 



On the other hand, from the virial relation 2T + U = 
between the halo kinetic energy T oc Ma 2 and its potential 
energy U oc GM 2 /r oc GM 5 ' 3 , the velocity dispersion takes 
the form 

a P ocMp /3 (l + z) 1/2 . (A4) 

Finally, the two preceding proportionalities can be cali- 
brated from the results of Wilson et al. (2001) on galaxy- 
galaxy lensing which give the following parameters for a halo 
harbouring an L* elliptical 

M, = 1.3 x 10 12 h' 1 M Q 
r, = 100 (50 to 200) h' 1 kpc 



a, = 255 (240 to 270) km s 1 . 



This leads to 



rp = r *U:J {1+z) 



(A5) 
(A6) 



r P ocMp /3 (l + z)- 



(A3) 



(Mp\ 1/3 n ,1/2 

With these approximations, the ratio between the for- 
mation rates -Rform inferred from the MPS model and the ki- 
netic theory (eqs. [3] and [Al], respectively) in the Einstein- 
de Sitter cosmology is constant and equal to unity for 
K — 1550 (see Fig. Al), except for large masses relative 
to the typical mass for collapse where the kinetic approach 
is no longer valid owing to the rarity of haloes. In this simple 
scale- free, n = power-law cosmology the behaviour of the 
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Figure Bl. Typical time interval between MMs vs. cosmic time 
for halo masses ranging from 10 11 to 10 15 M (top to bottom). 



MM rate as the product of two PS mass functions for mod- 
erate and small Mh can be seen from equation (2) by taking 
into account the corresponding explicit form of ct(Mh). Al- 
though this is not so obvious in the general case, the fact 
that one always finds very similar results to those plotted 
in Fig. Al allows extending that conclusion to any arbitrary 
cosmology. 



APPENDIX B: MBH AND HALO PROPERTIES 
Bl Growth of haloes 

As shown in Fig. (Bl), the time dependence of tmm(Mh, i) 
for a fixed mass Mh is very close to a power-law and the 
mass dependence for a fixed t is not far from this simple 
law, too. Thus, 



^(MH, t ) = r.(l)^ (^)" 



(Bl) 



with t*, a and b equal to 0.73, 0.71 and 0.15, respectively, 
M* = 10 11 M and t« = 1 Gyr. By bringing this general 
dependence into equation (4), we are led to 



dM H Ml ft, 
~~dT 



T- If J [MhW] " 



(B2) 



which easily integrates to 

M* 



M H (t) 



with B — (bt,)/(ar,). The set of functions {Mu_{t)} is a one- 
parameter family, with each member specified by its value 
M c at some reference time t c . 

All halo masses tend to zero as cosmic time tends to 
zero by virtue of equation (B3). On the other hand, when 



time tends to infinity, the second member of equation (B3) 
tends to 



' M, 
.mZ, 



M» 
Mc 



B 



(B4) 



The requirement that both members of equation (B3) are 
positive at all redshifts imposes a condition on M c , which 
cannot be arbitrarily large but has to be smaller than 



M 



top 



M„ 



B 



1/6 



2.7 x 10 1 



M Q . (B5) 



Otherwise, any halo of mass M c > M top at time t c would 
reach an infinite mass at some finite time. Conversely, given 
a fixed mass M c = Mh, the corresponding time of reference 
t c must be larger than the value t top satisfying 



Ctop — t* 



b ImT 



")' 



l/o 



(B6) 



B2 Growth of Massive Black Holes 

B2.1 Eddington-limited accretion regime 

According to the assumption (6), equation (5) writes 

dM BH Mbh [Mh jt),t] 
dt t mm [Mh (t) , t] 

and, taking into account equation (4), we get 



dM BH 



dM H 



Mbh Mh 

This leads, through equation (B2), to 



M BH (t) = Mbhc 



M u (t) 



M c 



(B7) 



(B8) 



(B9) 



and, from equation (6), we also get the solution AM B h(£). In 
equation (B9), Mbh c stands for M B h(M c , t c ), the boundary 
condition specified below. 

For any pair of physically acceptable values (Mm,ti) 
(i.e. with corresponding value of M c at t c , determined 
through equation [B3] for Mu(ti) = Mh;, satisfying the con- 
straint [B5]), the value M-Bn{M-m,U) is simply equal to the 
solution (B9) for t = U. 



B2.2 Supply-limited accretion regime 

Equations (5) and (9) and the time derivative of equation 
(8) lead to 



dM G _ M g [M H {t),t] 
dt eg Tmm [Mh (£),£] 

while equation (4) and the definition of M b leads to 

dM b M b [M H (t)] 



dt 



TMM[M H (t),t] 



(B10) 



(Bll) 



By denning the baryonic fraction in hot gas, 7(Mh,4) = 
M s (M H ,t)/M b (M H ,t), equations (B10) and (Bll) lead to 



TMM[M H (t),t] + 7 (t) W 1 , 



(B12) 
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with trivial peculiar solution 70 equal to (1 + e/eg) -1 ■ Then, 
by replacing equations (Bl) and (B3) into equation (B12) 
and solving the differential equation, one gets 



lit) = 70 + (7c - 7o) 



M H (t) 



1+-^ 

€ G 



(B13) 



with 7c the value of 7 at t c . 

As in B2.1, for any pair of physically acceptable values 
(Mm, U), the value j(Mjh, U) is given by the solution (B13) 
for t — U. The boundary condition 7 C = ~/(M c ,t c ) necessary 
to find this latter solution is actually unknown, as it requires 
to follow in detail the physics of baryons trapped in haloes 
of all masses since the dark ages. However, since we expect 
the change in accretion regimes from Eddington to supply- 
limited to take place when the baryonic fraction in hot gas 
becomes smaller than some threshold 7 C , we can assume that 
the regime changes at some time t c (or redshift z c ), for which 
all haloes have similar values of 7 equal to that threshold, 
i.e., y(M c ,t c ) = 7c- 

Therefore, for any given couple of values of z c and 7 C 
to be determined, one can infer the function ~/(Mu,t) and, 
from it and equations (7)-(5), all quantities of interest; first 

M s (Mn,t) =7(M H ,t)/ b M H (B14) 
M G (M H ,t) = [l--f(M H ,t)]f h M H , (B15) 

and then 

Mbh ( A/h , i) = £g Mg {Mh , i) (B16) 
AMbh (Mh , t) = e g M g (Mh ,t) . (B17) 

The values of these two latter functions at z c fix the bound- 
ary condition for the solution in the Eddington-limited 
regime, given by equations (B9) and (6). Let us finally note 
that by dividing equation (B17) by equation (B16) at z c , 
where the solutions Mbh(Mh,£) and AMbh(Mh,4) in both 
regimes match and, taking into account equations (B14) and 
(B15), we obtain 

■* = kh:- (B18) 

Thus, instead of z c and j c , we can use z c and k (for a given 
value of e g /eG) as the boundary conditions fixing the evolv- 
ing properties of MBHs. 



